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Abstract: We consider the problem of the limit of bio-inspired spatially extended neu- 
ronal networks including an infinite number of neuronal types (space locations), with 
space-dependent propagation delays modeling neural fields. The propagation of chaos 
property is proved in this setting under mild assumptions on the neuronal dynamics, 
valid for most models used in neuroscicnce, in a mesoscopic limit, the neural-field limit, 
in which we can resolve quite fine structure of the neuron's activity in space and where 
averaging effects occur. The mean-field equations obtained are of a new type: they take 
the form of well-posed infinite-dimensional delayed integro-differential equations with a 
nonlocal mean-field term and a singular spatio-temporal Brownian motion. We also show 
how these intricate equations can be used in practice to uncover mathematically the pre- 
cise mesoscopic dynamics of the neural field in a particular model where the mean-field 
equations exactly reduce to deterministic nonlinear delayed integro-differential equations. 
These results have several theoretical implications in neurosciences we review in the dis- 
cussion. 
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Introduction 

Brain's activity results of the complex interplay of different cells, in particular neurons, electri- 
cal cells that manifest highly complex nonlinear behaviors characterized by the intense presence 
of noise. Neurons form large population assemblies at the scale of which emerge reliable and 
adapted responses to stimuli. Such local neural populations, often termed cortical columns, 
have a diameter of about 50/im to 1mm, contain few thousands to one hundred thousand neu- 
rons and are in charge of specific functions [25] . The interaction of several columns at different 
spatial locations allow processing of the complex sensory or cortical information and support 
brain functions. Such groups of cortical columns organize on the surface of the cortex and form 
spatially extended structures called neural fields, the activity of which is precisely at the scale 
most usual imaging techniques (e.g. EEG/MEG, optical imaging) record relevant phenomena, 
and also correspond to anatomical information revealed experimentally. A paradigmatic exam- 
ple is given by the primary visual cortex of certain mammals. In such cortical areas, neurons 
organize into columns responding preferentially to specific orientations in visual stimuli, and 
display specific connection patterns [19, 3]. The communication between neurons is character- 
ized by a delay due to the transport of information through axons and to the typical time the 
synaptic machinery needs to transmit it. These delays have a clear role in shaping the neuronal 
activity, as established by different authors (see e.g. [7, 28]). In such structures, several highly 
populated columns interact, and the number of neurons in each column is orders of magnitude 
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higher than the number of columns (e.g., orientations) involved. A variety of important brain 
states rely on the coordinated behaviors of large neural assemblies and recently raised the in- 
terest of physiologists and computational neuroscientists. Among these, we shall cite the rapid 
complex answers to specific stimuli [32], decorrelated activity [11, 26], large scale oscillations [5], 
synchronization [20], and spatio-temporal pattern formation [14, 6]. 

The mathematical and computational analysis of the dynamics of neural fields relies almost 
exclusively on the use of heuristic models since the seminal work of Wilson, Cowan and Amari 
[1, 36]. These approach implicitly consider that averaging effects counterbalance the prominent 
noisy aspect of in vivo firing observed experimentally, and describe the mesoscopic cortical ac- 
tivity through a deterministic, scalar variable whose dynamics is given by integro-differential 
equations. This model was widely studied analytically and numerically, and successfully ac- 
counted for hallucination patterns, binocular rivalry and synchronization [22, 13]. Justifying 
these models starting from biologically realistic settings has since then been a great endeavor [4] . 

In this manuscript we undertake a rigorous analysis of neural fields. From the biological 
viewpoint, these are spatially extended cortical structures made of several highly populated 
neuronal ensembles (the neural populations) in charge of specific functions. From the math- 
ematical viewpoint, neural fields are adequately described as the limit of a set of interacting 
nonlinear stochastic processes (generally governing the neuron's electrical potential and related 
variables) gathering into different homogeneous populations at specific locations on the cortex. 
Neurons in each population have similar dynamics and communicate with neurons of different 
populations depending on the respective positions of the populations on the cortex and after 
a specific time delay. In what we will call the neural field limit, both the number of neurons 
and the number of populations tend to infinity so that the populations completely cover a 
continuous space (a piece of cortex or a functional space). 

This problem is evocative of statistical fluid mechanics and more generally interacting particle 
systems, and as such has been widely studied in mathematics and physics chiefly motivated by 
thermodynamics or fluid dynamics questions. In particular, the probability distribution of a 
typical set of particles in the limit where the total number of particles goes to infinity, and 
fluctuations around this limit where characterized for a number of models [24, 10, 31, 30, 29]. 
It was shown in several contexts that when considering that all particles have independent 
identically distributed initial conditions (chaotic initial conditions), then in the limit where the 
number of particles tends to infinity, the behavior of a few particles remains independent as 
time goes by, and all particles have the same probability distribution, which is the solution of a 
nonlinear Markov equation, often referred to as the McKean- Vlasov equation. The underlying 
biological problem motivates the introduction of a notion of spatial labeling of the (fixed) 
neurons, involving two mathematical aspects that were not covered in the literature. First 
is the fact that this induces the presence of infinitely many types of neurons (corresponding 
to the column neurons belong to), and second is the fact that since neurons communicate 
through the emission of electrical impulses transported at finite speed through the axons, space- 
dependent delays occur in the communication between two cells. These two aspects necessitate 
the development of the propagation of chaos theory towards infinite-dimensional functional 
settings that we aim at achieving in the present manuscript. We will show that in the neural 
field limit, the propagation of chaos property holds. Moreover, the activity is shown to converge 
in a certain sense towards the solution of a new object, a delayed integro-differential mean- 
field equation with space-dependent delays. This object has substantial differences with the 
usual McKean- Vlasov limits: beyond the presence of delays, the neural field limit regime is 
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at a mesoscopic scale where averaging effects locally to occur, but is fine enough to resolve 
brain's structure and its activity, resulting in the presence of an integral term over space. The 
speed of convergence towards the mean-field equations is quantified and involves two terms, 
one governing the averaging effect in each population and the second corresponding to the 
continuum limit. In the neural field regime, the limit equations are very singular, in particular 
trajectories are not measurable with respect to the space. These limits are very hard to analyze 
at this level of generality, and even numerical simulations are intricate in finite-populations 
settings [2]. However, in the type of models usually considered in the study of neural fields, 
namely the firing-rate model, we show in a companion article [33] that the behavior can be 
rigorously and exactly reduced to a system of deterministic integro-differential equations that 
are compatible with the usual Wilson and Cowan system in the zero noise limit. Noise intervenes 
in these equations a nonlinear fashion, fundamentally shaping in the macroscopic dynamics. 

The paper is organized as follows. We start in section 1 by describing the mathematical setting 
of the study, abstracting classical relevant neuronal models that are specified and reviewed in 
appendix A, and more general models are considered in appendix B. We then analyze the 
intcgro-diffcrcntial delayed McKean-Vlasov equations that will constitute our limit neural field 
equation in section 2 and demonstrate in particular their well-posedness, before addressing in 
section 3 the propagation of chaos property and convergence of the network equations towards 
the solutions of the mean-field equation. In section 4 we illustrate how this approach can be 
used in practice to analyze the effect of the parameters on the dynamics of the system in a 
particular example, reviewing some results of [33] afresh on a new example where noise, delays 
and spatial structure interact to shape the mesoscopic response of the neural field. The results 
of the mathematical analysis are then confronted to different recent experimental observations 
on collective dynamics of neural fields in the brain, and a few open problems of interest are 
discussed in the conclusion section 5. 



1. Mathematical Setting 

In all the manuscript, we are working in a complete probability space (fi, J 7 , P) endowed with 
a filtration (J : t) t satisfying the usual conditions. We consider a spatially extended network 
composed of N neurons, each neuron belonging to one of P(N) populations characterized 
by their locations (?!,••• ,rp(jv)) € r p ^ Ar ) on the cortex (or the feature space) T, a finite- 
dimensional compact set 1 . The state of each neuron i in the network is described by a d- 
dimcnsional variable X % e E := R d , typically corresponding to the membrane potential of the 
neuron and possibly additional variables such as those related to ionic concentrations and gated 
channels described in appendix A, and satisfy the network equations: 

1 P{N) 1 

1 ' 7=1 p(i)=7 7 

+ a(r)dWi (1) 

iWhen considering T as the cortex, it will be a compact subset of R 9 , q = 2 or 3, and when considering 
that populations are defined by the neuron's function, the shape of T can take different forms depending on the 
geometry of the feature space, for instance in the case of the primary visual area, neurons code for a preferred 
orientation of a visual stimuli that can be represented in the torus r = $ . 
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where f(r, t,x) :FxEx£h>£ governs the intrinsic dynamics of each cell, (W/) is a sequence 
of m-dimensional Brownian motions modeling the external noise and a(r) : T H> R rfxm a 
bounded and measurable function of r G T modeling the level of noise at each space location, 
and b(r,r' ,x,y) : T 2 x E 2 H> E the interaction function of a neuron located at r' with voltage 
y on a neuron at location r with voltage x. The function r(r, r') : T 2 i— > R + is the interaction 
delay between neurons located at r and those at r' which is assumed to be a regular function of 
its two variables. We assume that all delays are bounded by a finite quantity r. The quantity 
r a is called the location of the population and a is the population label. For a neuron i in 
the network, the population function p:Nn>N associates to a neuron i the population a it 
belongs to. The number of neurons in each population in a network of size N defines a sequence 
of population size {Ni(N), N P(N) (N)) (we hence have N 7 (N) = N) corresponding 

to the number of neurons in population 7 when the network size is equal to N. The number 
of populations P(N) and the number of neurons in each of these populations is assumed to be 
deterministic 2 . The interaction term presents a scaling factor P ^^ N ensuring the boundedness 
of the input received by neurons from population 7 to the other populations, a biological fact 
related to the brain function and to the finiteness of the resources available for the synaptic 
transmission. 

The different locations r~ l of the populations are related to the organization of the neurons on 
the space T. These locations are distributed according to a specific probability measure A on T 3 . 
The locations of the P(N) populations, (ri, • • • , r P rm) £ T P ^ N \ are assumed to be randomly 
and independently drawn in T according to the probability A(dr) in a different probability space 
(Q 1 , J 7 ', TP'). We will denote by £ the expectation over the realizations of the space locations 
(r a )- 




r 



Fig 1. A typical architecture of neural fields: cylinders represent neural populations as cortical columns spanning 
across the cortex. Neuron i (red population atri £ Y) receives a spike from neuron j (green population at rj G T) 
after a delay r(ri,rj) creating a current b(ri,rj, X%, X%). 

It is clear that the larger the number of populations, the smaller the mean number of neurons 
per populations. The number of populations will hence compete with the typical number of 
neurons per population and hence with averaging effects. In the present article, motivated by 
the fact that the number of neurons in each population is orders of magnitude larger than the 
number of populations (see e.g. [16]), we will make the following assumption, referred to as the 

2 It is easy to generalize to random population number and population size. 
3 In the example of the visual area VI, A is the uniform measure on S 1 
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1 P[N) 1 

^-PME^V)^ (2) 



In the case of an infinite number of populations, this assumption ensures heuristically most 
populations are made of a diverging number of neurons 4 . 

The parameters of the system are assumed to satisfy the following assumptions: 

(HI). f(r,t, ■) is uniformly Kf Lipschitz-continuous, 
(H2). b(r,r', •,•) is uniformly L-Lipschitz-continuous 
(H3). There exists a K > such that 

\b{r,r',x,z)\ 2 < K(l + \x\ 2 ) 

(H4). The drift satisfies uniformly in space (r) and time (i), the inequality: 

|/M,x)| 2 <c(i + M 2 ) 

(H5). The drift, delay, diffusion and coupling functions are regular with respect to space variables 
(r, r') 6 T 2 (at least measurable, in practice generally assumed continuous). 

Let us first state the following proposition ensuring well-posedness of the network system under 
the assumptions of the section: 

Proposition 1. Let (X t °) tG [_ T0 ] a square integrable process with values in E N . Under the 
assumptions of the section, there exists a unique strong solution to the network equations (1) 
with initial condition X° , which is square integrable and defined for all times. 

The proof of this proposition is a direct application of Da Prato [8] as used by Mao [23] , and 
essentially uses the same arguments as those of the proof theorem 2. The interested reader is 
invited to follow the steps of the demonstration of that theorem to prove proposition 1. 

We are interested in the limit of such systems as the number of neurons N goes to infinity, 
under the neural field limit condition. 



REMARK. Let us briefly bring some results from the analysis of finite populations networks (i.e. the 
case where P(N) remains finite as N — > oo), which can be seen as a particular case of the current setting 
under the assumption that X is a sum of Dirac masses. In that case, the neural field regime considers 
that the number of neurons in each population tends to infinity. Standard theory allows proving that the 
network converge towards P coupled McKean-Vlasov equations: 

1 p 

dX t (r a ) = f(r,t,X t {r a ))dt + o(r a )dW t a + - J2 Ez [&(»"«, r 7 , X t ( 

7 = 1 

4 If all populations have approximately the same number of neurons, each JVj. (N) will be of the order N/P(N) 
and the condition (2) is satisfied when P(N) = o(N). The condition also ensures the size of most populations 
tend to infinity. Indeed, for instance if all but one population contains just 1 neuron, the last population contains 
N - P(N) neurons, and the sum is equal to 1 - l/P(N) + 1/(P(N)(N - P(N))) > 1 - 1/P(N) which will not 
tend to zero. 
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where (W?) are P independent Brownian motions. This model can be seen as a discrete approximation 
of the continuous neural field. When the asymptotic number of populations is infinite, corresponding 
heuristically to refining the spatial discretization (or increasing the number of neurons) one is likely 
to face two main difficulties: (i) the network equations will involve an infinite number of independent 
Brownian motions, one for each space location, and (ii) it will involve a limit, as P goes to infinity, of a 
sum of the mean-field interaction terms (it is rather a simultaneous limit under the scaling property (2) ). 
The infinite number of independent Brownian motions is not a technical artifact, but a fact related to 
the very nature of the problem: distinct neurons are driven by independent Brownian motions whatever 
their respective locations on the neural field Y , and no continuity or measurability is to be expected. 

In order to handle the first point, we introduce a particular object, the spatially chaotic 5 
Brownian motion on T, a two-parameter process (t, r) G R + x T H> W t (r) such that for any 
fixed r G r, the process t i— > Wt{r) is a d-dimensional standard Brownian motion, and for 
r 7^ r' in T, the processes W t (r) and W t (r') are independent. This process is relatively singular 
seen as a spatio-temporal process: in particular, it is not measurable with respect to the Borcl 
algebra B(T) of T. This object, defined as a collection of independent Brownian motions, clearly 
exists. More generally, in what follows, we will qualify a process Ct( r ) of spatially chaotic if the 
processes (t(r) and Ct( r ') are independent for any r ^ r' . 

We will show that the network equations (1) satisfies the propagation of chaos property in 
the limit where N goes to infinity under the neural field assumption, and that the state of 
the network converges towards a very particular McKean-Vlasov equation involving a spatially 
chaotic Brownian motion. In detail, for almost all realizations of the spatial locations (t" 7 ,7 G 
N) i.i.d. with law A, the asymptotic law of neurons located at r in the support of A will be 
measurable with respect to (T, B(T)) and converge towards the stochastic neural field mean-field 
equation with delays: 

dX t (r) = f(r,t,X t (r))dt + a(r)dW t (r) + J m 2 [b(r, r', X t (r), Z t _ T{ry) (r'))] dX(r') dt (3) 

where {Wt(r))t>o, r er is a spatially chaotic Brownian and the process (Z) is independent and 
has the same law as (X). In other words, we will show that the law of the solution X t (r), 
noted m(t,r)(dy), is measurable with respect to B(T), and that the mean-field equation can be 
expressed as the integro-differential McKean-Vlasov equation: 

dX t (r) = f(r,t,X t (r))dt+ f f b(r,r' ' , X t (r),y)m(t-T(r,r'),r')(dy) d\(r') dt+a(r)dW t (r). 

Jr je 

which will also be written, denoting £ r > is the expectation with respect to the distribution of 
the population locations over T with distribution A(-), 

dX t (r) = f(r,t,X t {r))dt + a{r)dW t (r) + £ r , [E^[6(r, r' , X t (r), Z t _ r{ry) (/))]] dt 
Let us eventually give the Fokker-Planck equation on the possible density p(t, r, y) of m(r, t) 

5 We use the term chaotic in the statistical physics sense as understood by Boltzmann's in his notion of 
molecular chaos ("Stofizahlansatz"). 
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with respect to Lebesgue's measure: 

d t p(t, r, x) = -V, | (f(r, t, x) + j b{r, r' , x, y)p{t - r(r, r'),r' , y) dX{r')^j p{t, x) 

+ \^ x [\a{r)\ 2 p{t,x)}. (4) 

The mean-field equations (3) are of a new type: they resemble McKean-Vlasov equations but 
involve delays, spatially chaotic Brownian motions and an 'integral over spatial locations'. This 
is hence a very unusual stochastic equation we need to thoroughly study in order to ensure 
that these make sense and are well-posed. The existence and uniqueness of solutions to these 
equations is addressed in section 2, and the proof of the propagation of chaos and convergence 
of the network equations towards the solutions of that equations is addressed in section 3. 
REMARK. Note that the setting considered here, though relatively general, can be further extended using 
locally Lip schitz- continuous drift and state- dependent diffusion functions as certain neuronal models 
require. Such refinements do not modify the results and the principles of the proofs, but induce an 
important increase of complexity in the presentation, and are commented in appendix B. 

2. Analysis of the mean-field equation 

The mean-field equation (3) involves two unusual terms: a stochastic integral involving spatially 
chaotic Brownian motions and an integrated McKean-Vlasov mean-field term. 

Let us start by discussing properties of stochastic integrals with respect to a spatially chaotic 
Brownian Brownian. Considering A t (r) a J-" t -progressively measurable process indexed by r G T 
such that for any r € T we have 

t 

E[|A s (r)| 2 ] ds <oc, (5) 

it is trivial to see that for any r G I\ the process N t (r) := L A s (r)dW s (r) is a well defined, 

square integrable martingale with quadratic variation J Q |A s (r)| 2 ds. 

The possible solutions (X t (r)) t , r of the mean-field equation have a law belonging to the set of 
probability measures on the continuous functions of [— t, T] with values in the set of mappings 
of r in E. It is important to note at this point that similarly to the spatially chaotic Brownian 
motion, the solutions of the mean field equations are not measurable in (T, B(T)) since the 
solution considered at different space locations r and r' in T, namely X t (r) and X t (r'), are 
independent. 

Though trajectories of spatially chaotic processes are non measurable, their probability dis- 
tribution, defining a set of measures parametrized by r G T, might be measurable. This is a 
necessary property to make sense of the mean-field equations. Handling this subtlety necessi- 
tates to thoroughly define the space in which we are working and where the mean-field equations 
are well-defined. We define Z the set of random variables whose law is measurable with respect 
to B(T) (the random variable itself is not assumed measurable with respect to B(T)). For a 
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random variable (Z(r)) rS r m 2, with law p(r, dx) measurable, we define with a slight abuse of 
notations the !L|(r) norm on T by defining, for f a (fi', J 7 ', TP') random variable with law A: 

\\Z\\l m = M[£ r [\Z(f)\ 2 ]} = f f x 2 p{r,dx)d\{r) (6) 

J r J e 



where £ denotes the expectation on Q! '. This clearly defines a norm on random variables indexed 



by r £ r, when identifying processes that are A ® P-a.s. equal. We denote !L|(r) the set of of 



random variables in Z such that 11^11^2 ^ < oo. 
Example: 

(i). The spatially chaotic Browman motion at fixed time t has, for all r £ V, the law of a standard 
Brownian motion. This law, independent ofr£T, is hence measurable with respect to (T,B(T)). 
Moreover, it belongs to IL^fT) and has a norm equal to t. 

(ii). Another example is given by the variable Zr(r) = J Q T A s (r)dW s (r) where A is a function 
of R + x r measurable with respect to the a-algebra S(R + ) <g> B(T) and satisfying the condi- 
tion L J r j A s (r)| 2 d\(r) ds. The thus defined variable is not measurable with respect to B(F), 
but belongs to Z since this variable is a centered Gaussian process with measurable variance 
\A s (r)\ 2 ds, hence the law of Z(r) is B{T) -measurable. Eventually, Z £ IL^fT) with 11^11^2 = 

J T J r \A a (r)\ 2 d\(r)ds. 

We extend this norm to processes with values in IL^(r). For (Zt(r) )te[u,v] a stochastic process 
indexed by r £ T such that the law of Z t (r) is measurable with respect to B(T), we say that it 
belongs to M := M 2 ([u, v), ~L\(T)) if we have: 



II M 



:= £ f (B 



sup \Z s (r)\< 

s£ [u,v] 



< OO 



and this quantity defines a norm on A4 2 ([u, v], lL|(r)) where are identified the processes that 

are A-a.s. and P-a.s. equal for all times. 

Example: 

(i). The spatially chaotic Brownian motion on [0,T] belongs to M 2 ([0, T], IL^(r)) and has a norm 
equal to T thanks to the classical property that the supremum of the Brownian motion has the 
law of the absolute value value of the Brownian motion. 

(ii). The process Z t (r) = J* A s (r)dW s (r) introduced above belongs to M := M 2 ([0, T], ILj(r)) and, 
thanks to Burkholder- Davis- Gundy inequality, has a norm < 4 f r \ A s (r)\ 2 dX(r)ds. 

Now that these norms are introduced, we are in position to show the well-posedness of the 
mean-field equations: 

Theorem 2. For any ((t( r )> t € [ — r ;0]> r £ T) £ M 2 ([-t, 0], L|(r)) a square-integrable 
process, the mean-field equation (3) with initial condition £° has a unique strong solution on 
[— r, T] for any T > 0. 



Proof. As usually done for this kind of properties, we reduce the problem to the existence and 



uniqueness of a fixed point of a map $ acting on stochastic processes X in A4 2 ([— r, T], IL^(r)) 
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defined by: 



CoV) + Jo f(r, s, X s (r))ds + /* <7(r)dW a (r) 

+ /„* / r E z [6(r, r', X s (r), Z s _ T{ry) {r'))]d\{r') ds, 
C°(r) , tGt-r.O] 
(Z t ) =(X t ) independent of (X t ) and {W t {-)) 



t > 



The first question we may ask is whether if this function is well dchncd, and if $(X) defines a 
process that belongs to Ai 2 ([— t, T], L^fT)). The initial condition, Stieljes integral and stochas- 
tic integral with spatially chaotic Brownian motions are well defined as already seen. Moreover, 
under the assumption that X t (r) has a law measurable with respect to B(T), each of these terms 
have a law measurable with respect to B(T). The mean-field term is slightly more complex. Let 
us denote by m{t,r,dy) the distribution of X t {r). This term can be written as: 

b(r, r', X s (r),y)m(s — r(r, r'), r', dy)d\{r') ds 

r Je 

This quantity is well defined and has a distribution B(T) -measurable with respect to r since 
we assumed ra(i,r, dy) is measurable with respect to r as well as (r, r') H ► b(r,r' ,x,y) and 
(r, r') i — ^ r(r,r') (assumption (H5)). Let us now show square integrability of &(X). We have: 



<4 ||Co°| 



t 2 A (r) ^ 

T_ftT / (1 
Jo 



TC f (1 
./o 



\X S 



li? (r) 



|cr(r)| 2 ciA(r)ds (7) 



which is finite since XeM^. 

We may hence iterate the map $. We fix X a process in Ai 2 ([— t, T], L| (r)) and build 
the sequence X fc by induction through the recursion relationship X k+1 = &(X k ). We aim 
at showing that these processes constitute a Cauchy sequence in M. 2 ([— t, T], 1L| (r)), and we 
introduce the norm of the process up to time t: 



\X\ 



M, 



:= E 



£r sup \X.(r)\' 
\se\-T,t] 



= E 



sup \X s (r)\ 2 d\(r) 
r se[-r,i] 



We denote by (Z fe ) a sequence of processes independent of the collection of processes (X k ) and 



having the same law. We study the norm — X h 



Mi 



We decompose this difference into 
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the sum of three elementary terms as follows, for t £ [0, T] and r £ T: 
X k +\r) - X k (r) = I* {(f(r,s,X*(r)) - f(r, a, X^r)))} ds 

'^ z [b(r,r',X^r),Zl T ^ rl) (r'))] 

~ m z [b(r, r', X^(r), Z k s _ T{ry) (,•'))]) }dX(r') ds 

[(M z [b(ry,X*-\r),Zl T{ry) (r'))} 

- M z [b(r, r', X^(r), Z^ (r rl) (r'))]) }dX(r') ds 
=: A t (r) + B t (r) + C t (r) 
We hence obviously have: 

M* := \\X k +i - X k f M2t < 3 (\\Af Mi + \\B\? M * + ||C||^) 
We treat each term separately. We have: 
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E 



(Cauchy-Schwartz) < Tiff E 



( sup | / f(r,u,Xl{r))-f{r,u,Xl-\r))du\ 2 )d\{r) 
r v se[o,t] Jo 

\X k s {r)-X k -\r)fds)d\{r) 



< TK 



r x Jo 



which directly implies ||^4||Jt(2 < T K'j f Q M k 1 ds. The terms B t and Ct can be controlled using 
the same techniques. Let us for instance treat the case of C t . We have: 



\\C\ 



M 2 



E 



sup 

r se[Q,t] 



(Ez^ry,^ 1 ^),^^/))] 



-M z [b(ry,XtHr),Z^l {ry) (r'))})d u dX(r') 



dX{r) 



E 



{CS) < t 

It 2 J" 
(H2) < t I? 

T 2 Jo 
t 

<t L z I E 



E 



z[|6(r,r^X^ 1 (r),Z^ T(rV) (rO)-6(r,r^X^ 1 (r),Z^ (^ , J ,, ) (/))]| ) 
Efl^W-^-^r)! 2 ] dsdX(r)dX(r>) 



du 



dX(r)dX(r') 



(i 



ds =t L 2 I M^ds 



The term B t is treated exactly in the same manner and yields the inequality: 



\\B\\ 2 M2t <tL 2 [ Mt'ds 
Jo 
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All together we obtain, using the fact that for all k > 1, t £ [— r, 0] and r € Y we have 
X*+ 1 (r)=^t*(0 = C?(r): 



A"' / \/ r/.s (S) 

y + 2A 2 ), readily in 

M k < y y . . . y M o dai _ _ _ dSk < k_j_ m o (9) 

and M t ° is finite since we assumed X € A4 2 " and showed that §>(X) € A4 2 -. Routine methods 
starting from inequality (9) using Benayme-Tchebychev inequality and Borel-Cantelli lemma 
allow proving existence and uniqueness of fixed point for $ (see e.g. [27, pp. 376-377]), and 
that this fixed point is adapted and almost surely continuous. Moreover, this process has a 
measurable law with respect to (Y,B(Y)) as a limit of measurable laws, and since it satisfies 
$>(X) = X, it is a solution to equation (3). Let us eventually show that it belongs to M.\. We 
stop the process at time 

r„ = inf{f > 0, 11X11^2 > n}. 
Using inequality (7) and the fact that X — $(X) we have: 



with K' = 3T(AT 2 + 2A 2 ), readily implying: 



(llCo L i(r) + t (c + k + 4 |M| Li(r) ) + t(c + k) j A " \\x\\ Ml ds^j 



ensuring by Gronwall's lemma that r n a.s. tends to infinity as n — > oo, and eventually, letting 
n — > oo, that the solution has a finite norm in 

Proving uniqueness of the solution using equation (8) is then folklore. 

□ 

Now that we proved strong existence and uniqueness of solutions for the mean-field equations, 
we now turn to showing that network equations indeed convergence in law towards this solution, 
and that the propagation of chaos occurs. 



3. Limit in law and propagation of chaos 



We are now in a position to prove the main result of the manuscript, namely the convergence 
in law of the solutions of the network equations (1) towards the equations (3), and the fact 
that the propagation of chaos property occurs. To this end, we consider that the network 
equations have chaotic initial conditions. In details, let (Ct( r )) G Al 2 ([— r, 0], IL 2 (r)) a spatially 
chaotic stochastic process, i.e. a stochastic process such that for any r/r', the process ((t(r)) 
is independent of (C t °(V')). We consider that the initial condition of different neurons in the 
network are independent, and the initial condition (Q) g A^ 2 ([— r, 0], 1L 2 (Y)) for neuron i in 
population a, is equal to (Ct( r a)) € M 2 (C T ). 

The classical coupling argument cannot be directly applied here. Indeed, the usual argument 
is based on the fact that we are able to define the solution of the mean-field equation through 
the use of the same Brownian motion and with the same initial condition as one of the neurons 
(or particles). This is no more the case because individual neurons are governed by finite- 
dimensional Brownian motions and the mean-field equation by a spatially chaotic Brownian 
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motion. Notwithstanding, an argument based on a slightly more subtle couplings holds. In 
details, let us consider neuron i G N of the network, in population a at location r a 6 V. 
Denote by {WD the Brownian motion governing the evolution of neuron i in the network and 
C € A4(C T ) the initial condition of the network. We aim at defining a spatially chaotic Brownian 
motion Wl(r) on R mxd such that the standard Brownian motion (W| (r a )) is equal to (W^), and 
proceed as follows. Let (W t (r)) tS [ ,T],rer be a m x d-dimensional spatially chaotic Brownian 
motions independent of the processes (W*). The processes: 

(Wi(r)) = (W t (r)) r^r a 
{Wl{r a )) = (Wi) 

is clearly a spatially chaotic Brownian motions, and will be used to construct a particular solu- 
tion of the mean-field equations. In order to completely define a solution of the mean-field equa- 
tions, we need to specify an initial condition, and aim at coupling it to the initial condition of 
neuron i. To this end, we define a spatially chaotic process {(t{r)) G M 2 ([— r, 0], IL|(r)) equal in 
law to (Ct°( r )) and independent of Q, and define a coupled process ((I'' (r)) G Ai 2 ([— r, 0], JL\(T)) 
as: 

ftf°(r)=C?(r) ^r Q 

ld'°(r Q ) = a- 

Here again, it is clear that this process is spatially chaotic, i.e. that for any r ^ r', the processes 
Ct°(r) and Ct'°( r ') are independent, and that Ct'°( r ) has the law of Ct°( r )- 

Now that these processes have been constructed, we are in a position to define the process 
(XI) as the unique solution of the mean-field equation (3), driven by the spatially chaotic 
Brownian motion (W'l(r)) and with the spatially chaotic initial condition (C t '°{r))\ 

( dX\ (r) = /(r, t, Xl(r)) dt + J r M z [b(r, r', X*(r), Z t _ T{r ^ (r'))\d\(r') dt 

+a(r)dWi(r) , t>0 
( Xl{r) =C°(r) , te[-r,0] 

k (Z t ) =(X t l ) G X independent of (Xj), {W l t {-)) and •)) 

The same procedure applied for all j G N allows building a collection of independent stochastic 
processes (X^ (r))j =1 jy G A / J 2 ([— r, T], 1L 2 (r)) such that all neurons j in population a have the 
same law as (X(r a )). Let us denote by m(t,r) the probability distribution of X t (r) solution 
of the mean-field equation (3). As previously, the process (Z t (r)) generically denotes a process 
belonging to M 2 ([— t, T], Jj\(T)) and distributed as m. 

Let us fix I G N* and • • • ,ii) a collection of neuron indexes such respectively belonging 
to populations located at (n, • • ■ , r^) (possibly identical). We now prove the almost sure con- 
vergence of a collection of processes (Xl k '" ,k = l...l) towards (X\ k (r k ), k = 1 . . . I), implying 
its convergence of the law towards the chaotic distribution m(t,r\) ® ... ® m(t, r k ) as ./V goes 
to infinity. We start by proving this property for / = 1 before extending that result to I > 1. 

Theorem 3. Let i G N a fixed neuron in population a. Under the assumptions (HI )-(H5) and 
the neural field assumption (2), for almost all realizations of the population locations (r a ,a G 
IN), the process (X t 1 ' ,t <T) solution of the network equations (1) converges in law towards the 
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process (X t (r a ),t < T) solution of the mean-field equations (3) with initial condition (Ct( r )) 
and moreoveover, the speed of convergence is given by: 



£ E 



sup \Xt' N - Xl(r a )\ 2 )=0[<b(N) + -±-\ (10) 



-T<S<T 



P(N) 



REMARK. We recall that£ denotes the expectation on the distribution of the space locations (rk)k=l'--P(N) 
and<B(iV) = ^E^ ) K^)- 

Proof. The proof is based on evaluating the distance E[sup_ r<s<T \X' l s ' N — X l s \ 2 }, and breaking 
it into a few elementary, easily controllable terms. A substantial difference with usual mean-field 
proofs is that we need to prove a convergence in the infinite-dimensional space 1L\(T), and that 
the interaction term in networks equations consists of a sum over a finite number of populations, 
whereas the effective interaction term arising in the mean-field equation is an integral over T. 

In all the demonstration, we will generically denote rp £ T the location of population f3 £ 
{1, • • • , P(N)}. We use the following elementary decomposition (each line of the righthand side 
corresponds to one term of the decomposition, A t (N) — E t (N)): 



Xi-Xi(r a )= f (f(r a ,s,Xl)-f(r a ,s,Xi(r a )))ds 
Jo 

\ I 7 =1 Jo T j=l 
i P(N) t i «i 

+ M E / T^E^'^^),^-.^))-^^ 



rf.s 



1 /■ / 1 - - - \ 

+ ^7M E / (^E fe ( r «' r 7'^( r «)'^-r(r cv ,r 7 )( r 7)) ~ [6(r Q , r 7 , X](r Q ), Z s _ T(r<i (r 7 ))] J 

\ I 7=1 "'O T j=l 

1 PW /■*/ /■ 

+ P(A7) E y o (Ez[&(r Q ,r 7 ,A:(r Q ),Z s _ T(ra ,^^ 



7=1 J ° 

=: A (AO + B t (N) + C t (N) + D t (N) + E t (N) 

Due to the exchangeability of neurons belonging to the same population, the probability dis- 
tribution of these terms do not depend on the particular neuron i considered but only on the 
population it belongs to. The terms A t (N), B t (N) and C t (N) involve the Lipschitz continuity 
of the functions involved, the term D t (N) correspond to averaging effects (mean-field limit) at 
single populations levels and the term E t (N) corresponds to the continuous limit. The terms 
A t (N) through Ct(N) are treated using the Lipschitz continuity of the functions involved. Using 
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Cauchy-Schwarz (CS) inequalities, we easily obtain: 



0<s<t 



E[sup \A s {N)\']<KiT I E[ sup \X^ N - Xl[r a )\']ds 



— T<U<S 



E[sup |S s (iV)| 2 ] < TL 2 E[ sup |^ -X;(r a )| 2 ]d S 



— T<U<S 

t 



E[ sup |C s (iV)| ] < TL / max E[ sup - X£(r pW )|^]d S 



0<s<t 

Let us for instance treat the case of B t (N) 

P{N 



3=1— N - T <u<s 



E[ sup |S s (iV)| 2 ] = 



1 



0<s<t 



P(N} 



-E 



, P(A,) r s i ^ / 

0<a<tl ^ Jo ^7^ «-r(ra,r,) 



(CS)< 



T 



'7 



du\ 



p(n) ^ y iv 7 ^ 



E 



6(r„, r 7 , Xl< N , X'i^ ^) - 6(r a , r 7> JC„, X J _ r(r 



((H2)) < TL- 



< TL Z 



Ml 



X 



,N 



XI ds 



E 



sup \X^-Xl 

—T<U<S 



ds. 



The mean- field term D t (N) involve the difference between an empirical mean of a function of 
processes and an expectation term, and all have bounded second moment thanks to theorem 2 
and assumption (H3). We have using (CS) inequality: 

P(N) t N 1 

E[sup \D S (N)\ 2 } < —— £ / E[— ^6(r a ,r 7 ,Xi,xf_ T(rairT) )-E z [6(r Q ,r 7) ^,^_ T(i 

0<s<t v / 7= i "^0 7 j — l 

and hence involves an expectation of type: 

N, 



E[ 



i Nl 1 ^ r 

\ — Y,Q(xi,xD-®zMx l s ,zm 2 } = jpJ2 ^[(&(xi,xi)-M z [e(xi,z])}) 

7 j=l T fc,Z=l 



iV 2 

(e(x:,x 5 fe )-E Z [e(x:,z 7 )]) 



where 0(x,y) — b(r a ,r^,x,y). Routine methods allow to show that all the terms of the sum 
corresponding to indexes j and k such that the three conditions j : ^ i, k ^ i and j ^ k are 
satisfied are null. One simple way to show this property consists in writing the expectations as 
integrals with respect to the measure m(t, r a ) and observing that all terms annihilate. Therefore, 
there are no more than 3 AL non-null terms in the sum (in the case a = 7 there are just JV 7 
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non-null terms), and moreover, all of these terms are uniformly bounded. The terms related to 
indexes j = k =/= i satisfy the inequality: 



E 



\e(xi,xi)-M z [e(xi,z2)] 



< 2E 



\e(xi,xi] 



\M z [Q(Xi,Z] 



< 2 + E[|X*| 2 ]) + E 

< 4K(1 + C'(s)) 



M z [JK(l + \Xl\i)] 



with C'(s) given by theorem 2. The terms related to the cases j = i (or symmetrically fc = i) 
are bounded by the same constant, since we have for all k such that p(k) = a, using Cauchy- 
Schwartz inequality. We note C = 4K(1 + C'{T)). We hence conclude that: 



E[sup \D S (N)\ 2 } <T 2 C-±— 

0<s<t ) 



P(N) 

E 

7=1 



3i\L 



1 



N 2 

7 



< 3T 2 C 



P(N) 



P(N) 

E 

7=1 



1 



3T 2 Cs(N) 



It hence only remains to control the term E t (N) corresponding to the difference between 
an integral over the space T weighted by the density dX(r) and a sum, weighted by 1/P(N) 
of the same integrand at P(N) discrete values (r 7 ) € T N independently drawn in T with the 
probability density d\(r). This sum hence resembles to a Monte-Carlo approximation of the 
integral term, and we now show that our sums over populations converges for almost all choice 
of (r 7 ) G T N towards the integral, using an argument similar to the one we just used to control 
D t (N). In details, we show that £ (E[sup 0<;j<t |-E S (./V)| 2 ]) converges towards 0, using the same 
method as used for the convergence of the mean-field term. Let us denote for the sake of 
compactness of notations F(s,r,r') the expectation Ez[d(r, r' , Xl(r), Z S _ T( > . r <)(V)]. 

We have: 



£(E[ sup \E S (N)\ 2 }) <T £ ( E 

0<s<t Jo 



1 



P(N) 



P(N) 

E 

7=1 



F(s,r Q ,r 7 ) - £ r , [F(s,r a ,r')] 



ds 



Similarly to what was done for the term D t (N), since £ r i [F(s, r Q , r')] is precisely the expectation 
of F(s, r a , r 7 ) under the law of r 7 over which the sum is taken, developing the squared sum into a 
double sum over populations (say, 7 and 7'), it is easy to show that, because of the independence 
of the r 7 , that all terms that do not correspond to 7 = 7', 7 = a or 7' = a vanish, leaving 
less than 3P(N) possibly non-null terms, and these terms are uniformly bounded. Indeed, for 
r 7 = r-yi (the case r 7 = r a is treated in the same manner), we have: 



£(E[\F{s,r a)ri ) - £ r >[F(i 



2 }) < 2£ (E[|F(s, r a , r 7 )) | 2 + \£ r > [F(s, r a , r'J] | 2 ]) 



<2£(E[|F( S ,r Q ,r 7 )| 2 
< 4K(1 + C(s)) 



■£ r ,[\F(s,r a ,r')\ 2 }}) 



implying eventually that: 



f(E[sup \E S (N)\ 2 ])< 4 j^r(l 

0<s<t P{J-V ) 



C{T)). 
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All together, we hence have: 



E[ sup \Xl' N - Xg(r 7 )| 2 ]) < K' 

Q<s<tATu 



valid for all i G M, hence we have: 



£[max E[sup \X* 3 ' N - Xj(r 7 )| 2 ])] 



tATu 



max E[ sup \X^ N - Xi(r pU) )\ 2 ]ds 

3=1-N _ T < n < s 

+ Cice(7V)+E[ sup \E S (N)\ 2 
a<s<t 



-1—N 



0<s<t 



< K' 



tAT V 



max . E[ sup \Xi> N - Xi(r p(j) )\ 2 } 

3=1— N _ T < n < s 



ds + d(E(iV) + 



Co 



P(N) 



where if' = AT (K 2 f + 2L 2 ), C x = 12 T 2 C and C* 2 = 16T 2 J ftT(l + C(T)) neither depend upon N 
nor in the particular neuron considered. By Gronwall's inequality, we hence obtain: 



max E 

j=l— N - T <S<tAT,j 



sup \Xi N -Xi(r,)\ 2 )) 



< (Ci<b(N) + 



c 2 

P(N)J K 1 



which ends the proof. 



□ 



Corollary 4. Let I € IN* and /ii Z neurons ■ ■ ■ ,ii) S N*. Under the assumptions of theo- 



rem 3, the law of (X. 



il,N 



,x. 



ii,N 



for almost all realization of the population locations (r a , a € N). 
Proof. We have: 



< t < T) converges towards mt{r p u x )) ® ■ ■ ■®mt{r p (i l )) 



£ (e 


sup 




■ ,x^ N )-{xt\- 


■ ,x?) 


2 




-T<t<T 











< 



k = l 



^E 


sup 




-T<t<T 



x 



K'T 



h,N 



.X 



il,N 



,-t <t<T) 
, X\ l , — t < t < T) whose law is equal by definition to 

□ 



P{N), 

which tends to zero as N goes to infinity, hence the law of (X 
converges towards that of (X l t 
m (t, r p(n)) ® • • • <8> m(t, r p(il) ). 

IMPORTANT REMARK. The speed of convergence towards the mean-field equation is hence governed 
by <b(N) and 1/P(N). In the case of a finite number of populations, the speed of convergence is hence 
driven by the size of the smallest population. In the infinite population case, the speed of convergence 
towards the mean-field limit is a balance between the averaged number of neurons in each population 
through the term <b(N), and the total number of populations through the term 1/P(N). The first term 
quantifies the speed at which averaging effects occur in the network, and is related to the averaged 
inverse number of neurons in each population. The other term controls the convergence of the interaction 
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related to all populations towards an effective interaction term given by an integral over Y of mean-field 
interactions, i. e. convergence of finite-populations networks towards their continuous limit. For networks 
with homogeneous population sizes, <b(N) will be approximately equal to P(N)/N . The optimal network 
size ensuring the fastest convergence in that case hence corresponds to P{N) ~ vN (minimizing the 
functional x t-¥ P(x) /x + l/ P(x) ), and in that case the convergence will be ml/ vN, and we conjecture 
that this speed of convergence is optimal (though we did not achieve to prove it). This convergence speed 
is hence very slow compared to finite-size networks and usual mean-field limits in which the speed of 
convergence is of order 1/N . 

4. Neural Fields Equations in Action 

It is folklore to notice that McKean-Vlasov limits have dynamics very complex to analyze. Very 
refined methods are generally set up to analyze the behavior of the system in the mean-field 
limit, such as entropy methods or spectral methods (see e.g. [34]). This statement could be 
even more true in our spatialized context, and the present, general approach might appear to 
be bound to remain formal. 

Fortunately, for relevant neuroscicncc applications, it happens that solutions of these equa- 
tions are not out of reach. This is the topic of a companion article [33] where is considered 
networks of firing-rate neurons (see appendix A. 2), the neuronal model usually considered for 
neural fields analysis. Let us briefly review here the main results of that article and concretely 
use the proposed approach to analyze the dynamics of a simple network. 

Considering firing-rate neurons, we show in [33] that the solutions of the mean-field equa- 
tions are Gaussian processes when the initial condition also is (and equilibria are Gaussian), 
and that their mean M(r,t) and standard deviation v(r,t) (fully describing the process since 
the covariance is a simple function of these two quantities in that case) reduce to the set of 
deterministic delayed integro-differential equations: 

'd t M(r,t) =-^M(r,t) + I(r,t) 

+ jT r J(r, r')F(r', M(r', t - r(r, r')), v(r',t - r(r, r')))X(r')dr' (11) 

A< r ^) = -^M) + ^M 2 

where F(r, x, y) denote the expectation of S{r, U) for U a Gaussian random variable of mean x 
and variance y, and can be made explicit for particular choices of sigmoids S. These equations 
are consistent with the heuristically derived extremely widely used Wilson-Cowan models for 
finite-populations neural assemblies [35, 36] in the limit where noise levels vanish. These equa- 
tions are shown to be well-posed, and address the existence of spatially homogeneous solution 
in law (i.e. solutions whose law do not depend on r £ T). The analysis of the moment equations 
are shown to grant access to the dynamics of the network, and this is developed in details for 
some biologically motivated models. In that article, the choice of the parameters, driven by 
biological constraints, did not reveal any qualitative effect of the delays on the solutions except 
during transient phases. 

In order to illustrate how the use of the present approach can be used to uncover the dynamics 
of the neural field, we proceed to the analysis of single population network with inhibitory 
interactions (i.e. negative interactions), a case that was not treated in [33] and which will turn 
out show a particularly rich variety of behaviors as a function of delays. 
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To this end, let us fix the parameters of the system. We consider T = S 1 the 1-dimensional 
torus, and A the uniform distribution on it. We consider that S(r,x) = J^ x e~ x ' 2 /\^2tt =: 
erf(gx) , 6(r) — 1 and a independent of r, and can easily show by integration by parts that 
F{x,y) = evf(gx / y/ 1 + g 2 y)- We further fix J(r, r') — Je^ r ^ r ^ & (8 represents the typical 
connectivity length in the neural field) and r(r, r') = \r — r'\/c + r s (c represents the speed of 
transmission in the neural field and t s the typical transmission time of the synapse) . 

Since -F\0, y) = for any y £ R, the Gaussian solutions with zero mean and standard 
deviation c 2 /2 are stationary solutions of the system that are spatially homogeneous in law. 
Characterizing the stability of this solution consists in analyzing the characteristic roots equa- 
tion of the linearized system around the spatially homogeneous stationary solution. Computing 
the eigenvalues of the integral convolution operator similarly to [33, Section 3.1], we obtain the 
dispersion relationship: 

c 

for k € 7L and F' a = , 9 -r=. The spatially homogeneous equilibrium is stable if and 

V 1 +9 2 * J ^ 

only if all solutions £ to the dispersion relationship (characteristic roots) have negative real 
parts. A Turing bifurcation point is defined by the fact that there exists an integer k such that 
5i(£) = 0. It is said to be static if at this point 9(C) = 0, and dynamic if 3(0 = tok 0. In 
that latter case, the instability is called Turing-Hopf bifurcation, and generates a global pattern 
with wavenumber k moving coherently at speed tOk/k as a periodic wavetrain. 

Possible Turing Hopf bifurcations hence arise when there exists Uk > such that: 



icj k + 1 = F^e-^'Zkiu}) 



with Zk(uj) = ^\?v~?*o ^t) i which yields bifurcation curves (parametrized by cj) in the parameter 
space: 

2 / i , JVI^MI 



a 2 =^1-1 



5 2 V 2tt(1+w 2 ; 
— arctan(w) + Arg(FQZfe(cj)) + 2m7r 

U! 

This provides a curve of Turing Hopf bifurcations corresponding to transitions from stationary 
independent solutions to perfectly synchronized independent solutions, as displayed in figure 
Fig. 2 and 3. In Fig. 2 we display the bifurcation curve in the parameter space (cr, r s ) for a 
specific set of parameters. This curve has a convex shape. Small enough delays hence correspond 
to stationary solutions. Increasing delays yields periodic activity, which disappears as noise is 
increased. This example shows the importance of delays in the qualitative dynamics of the 
neural field. The typical connectivity length also shapes the qualitative dynamics of the neural 
field, as shown in Fig. 3. This variety of behaviors correspond to bifurcations corresponding to 
a wavenumber fc = 0, and correspond to spatially homogeneous solutions. Non-trivial spatial 
structures can be searched for considering non-spatially homogeneous initial conditions. In 
this case, a number of complex spatio-temporal behaviors can appear, such as the metastable 
polychronization shown in Fig. 3 where the neural field splits into two clusters oscillating in 
antiphase during very long transient periods before a sudden synchronization of the whole 
neural field. 
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(a) Bifurcation Diagram (b) cr = 0.1, t s = 0.4 (c) cr = 0.1, r a = 0.5 (d) a = 0.3, t s = 0.5 



Fig 2. Turing-Hopf bifurcations and delay-induced synchronization, J = —3, g = 3, S = c = 1. faj.' bifurcation 
diagram, shows a transition from stationary to periodic activity as delays are increased (b)—¥(c). When noise is 
increased, synchronization is lost (c)—> (d). (b)-(d): spatio-temporal dynamics as a function of space (abscissa) 
and time (ordinate) 




(a) 5 = 1 (b) <5 = 5 




(c) Antiphase (d) Antiphase (e) Synchronization 



Fig 3. Spatial effects: increasing 6 destroys the synchronization. For 8 = 1 (case (c) of Fig. 2), choosing non- 
spatially homogeneous initial conditions yields to complex situations, for instance an antiphase synchronization 
during long transients (bottom row), (c): t £ [0,200], (d): orange: M(t, 0.1), blue: M(t, 0.9), black: M(t, 0.5), 
(e): t £ [600,650]. The synchronization becomes visually perfect for times above 1500. 
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5. Discussion 

In this paper, we addressed the problem of the asymptotic behavior of networks composed of 
a large number of neuronal assemblies in a particular asymptotic regime, the neural-field limit. 
We took into account a number of specificities relevant to neuronal dynamics: intrinsic noise at 
the level of each neuron, the spatial structure and propagation delays. We demonstrated that 
for a relatively general class of models, that includes the most prominent models in neuroscience 
(reviewed in appendix A), the propagation of chaos property took place, and showed convergence 
of the mean- field equations towards mean- field equations of a new type, analogous of the classical 
McKean-Vlasov equations but including delayed interactions, spatial integration term and a 
singular spatio-temporal stochastic process: the spatially chaotic Brownian motion. 

The question of the scale at which relevant phenomena occur is essential to the modeler. 
Descriptions coarser than our neural field limit, for instance those involving finite numbers of 
populations, correspond to cases where our measure A is a sum of Dirac masses. This case 
can be seen as a particular case of the present analysis, and hence the propagation of chaos 
occurs and network equations converge towards mean-field equations that correspond to a finite 
system of delayed McKean-Vlasov equations. In contrast, scales finer than the neural field limit 
(taking for instance into account possible individual heterogeneities between neurons) are not 
covered by the analysis and seem relatively hard to understand. It is likely that the dynamics 
of such networks will be considerably distinct from that of networks in the neural field regime. 
The neural-field regime seems particularly well suited to describe the activity of large neuronal 
assemblies, since it was observed that population sizes are orders of magnitude larger than the 
total number of populations [16]. Moreover, it seems to be at the scale of biological recordings 
and phenomena such as the emergence of patterns of activity in the cortex. We illustrated how 
such an analysis could be rigorously developed in simple example in section 4. More relevant 
states may be analyzed with this model, since the usual heuristic equations that were successfully 
used in a number of situations [13] are compatible, in the zero noise limit, with our equations, 
and the rigorously derived model will shed new light on the role of noise in such neuronal 
systems, but also on the individual behaviors of neurons. For instance, the propagation of chaos 
property ensures that finite sets of neurons are independent in the neural field limit. This result 
contradicts the classical view considering that since neurons of the same population are highly 
connected and receive similar input, their activity shall be correlated. However, with recent 
experimental findings using high quality recordings [11, 26] showed that levels of correlations 
between two neurons (of the same population or not) were extremely small, way below what 
was usually considered. The propagation of chaos hence offers a universal explanation to this 
phenomenon. 

A number of open questions remains widely open in the theoretical understanding of the 
behavior of neural fields and large-scale neural networks. For instance, a particularly interesting 
phenomenon is the plasticity of neuronal connections. Considered constant and homogeneous in 
the present manuscript, it happens that the synaptic coefficients describing pairwise interactions 
between neurons evolve, very slowly, as a function of the network activity and in particular as 
a function of the correlations between the activity of pairs of neurons. This kind of phenomena 
was never considered in the mathematical literature, and seems relatively rich. In particular, 
this mechanism can break the propagation of chaos property and yield weakly correlated states. 
This is a problem we are currently investigating. 
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Appendix A: Neuron models 

For the sake of completeness we quickly review in this appendix different classical neuron mod- 
els motivating the present study. This appendix takes a mathematical viewpoint, is obviously 
very selective and lacunar. The interested reader will find more details in classical neuroscience 
textbooks, e.g. [21, 12]. Basically, neurons are electrically excitable cells whose activity, mea- 
sured through the voltage of the cell (difference of electrical potential between the intracellular 
and extracellular domains), is governed by ionic transfers through specific proteins (ion voltage- 
gated ion channels) located on the cellular membrane. We present here detailed neuron models 
(appendix A.l) that approximate the biophysics of ion channels, and firing-rate models (ap- 
pendix A. 2) that reproduce qualitatively the dynamics of the firing rate of neurons and that 
are used in the application section 4. 

A.l. Hodgkin- Huxley and Fitzhugh-Nagumo models 

Probably the most biologically relevant, versatile and precise neuron model is the Hodgkin- 
Huxlcy (HH) model [18]. This model describes the membrane potential v of a neurons as a 
function of the dynamics of several ionic currents that enter or exit the cells through voltage- 
gated channels. The mathematical description we choose here involves Langevin approximation 
of the random proportion open of ion channel (see e.g. [17] and references therein). The pro- 
portion of open channels satisfies in that model a stochastic differential equation: 



where Wf are independent standard Brownian motions, A x {v) and B x (v) are smooth bounded 
functions accounting respectively for the opening and closing probability intensity of a given 
channel and x( x ) is a function vanishing outside [0, 1] to ensure that the variables x remain in 
[0, 1] (since these variables describe proportions). Generally, three ionic currents (and channels) 
are considered: potassium (to), calcium activation (n) and inactivation (h), and Ohmic leak 
current, II (carried by Cl~ ions). Considering that the neuron receives an external current 
composed of a deterministic part I(t) and a white noise with standard deviation a ex t, the 
voltage is governed by the equation: 



This model satisfies assumptions (HI) and (H4) used in the general theory, since though poly- 
nomial nonlinearities arise in the dynamics, the boundedness of the variables (n, to, h) ensure 
Lipschitz continuity and linear growth. The assumption (H5) is not satisfied since the noise 
depends on the state of the neuron. This refinement does not make the proofs substantially 
more intricate as discussed in appendix B. 

The HH model is often too complex for practical purposes and several reductions were 
proposed. A particularly interesting on is the Fitzhugh-Nagumo (FN) bidimensional model [15] 
capturing from the biological viewpoint the most prominent behaviors of the Hodgkin-Huxley 
model. From the mathematical viewpoint, it is important to specify this model since that 



dx t = (A x (v)(l -x)- B x (v) x) dt + ^A x (v)(l - x) + B x (v) x X (x) dWf 




(12) 
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model does not satisfy assumptions (HI) and (H4), and motivates the additional mathematical 
developments of appendix B. This model describes the evolution of the membrane potential 
variable v and a slower recovery variable w, through the equations: 

d,V t = (P( Vt ) -W t +I)dt + G V dW? 

dw t = a (b vt — Wt) dt + a w dW™ 

where P(v) — v(l — v)(v — o), generally chosen f(v) — v — v 3 . 

The state of the neuron X in our abstract model (1) in the HH model is given by (v, n, m, h) 
and for the FN model by (u, w) and their intrinsic dynamics is enclosed in the functions / and 
9- 

The communication between neurons is maintained by two possible types of synapses: elec- 
trical or chemical. Electrical synapses, in charge of rapid and stereotype signal transmission, is 
done through direct contact of the intracellular domain of the two communicating cells through 
specialized protein structures called gap-junctions. The ions passively flow from one neuron to 
the other: the interaction is not delayed and the current produced by neuron j on neuron i is 
equal to Jij(vj — v\) where Jij is called the synaptic conductance (this defines our interaction 
function b in the abstract model (1)). When including the dependence on v\ in the drift function, 
the interaction function J^j Jij v t clearly satisfies assumption (H2) and (H3), and (H5) as soon 
as the dependence of Jy with respect to space is sufficiently regular. The chemical synapse is 
the most common type of interconnection. When a spike is fired from a pre-synaptic neuron j, 
it is transported through the axons to the synaptic button where it is transmitted to neuron i 
through a complex process of release of neurotransmitter (from j) binding to specific receptors 
on neuron i. The transmission takes a time r.y of the order of a few millisecond. Similarly to HH 
ion channels dynamics, the proportion of open neurotransmitter channels y l has the dynamics 
(see [9]): 

dy{ = (AS{v*)(l - y j (t)) -Dyi(t)) dt + oytfrf) dW(' v . 

with S is a smooth sigmoidal function. In our abstract model, the variable y % is added to the 
state X 1 of neuron i and the functions / and g take into account that dynamics. The synaptic 
current induced at time t on neuron i by the arrival of a spike from neuron j (fired at time 
t—Tij) is equal to Jijy J (t—Tij)(v l (t) — v rcv ) governing our interaction function b clearly satisfying 
assumptions (H2) and (H3), and (H5) as soon as the dependence of Jy with respect to space 
is sufficiently regular. 

The synaptic efficacies Jij of electrical or chemical synapses are given by the connectivity of 
the cells. Such functions are generally considered continuous functions J{ri,rj) depending on 
the population of i and j. 

Putting all these elements together and assuming that all the parameters of the equations 
only depend on the neural populations of the cells involved, we can write the equation of a 
network of FN neurons with chemical synapses, external and synaptic noise: 

'dvl = (P(vl) + P(t) + E^i.i^ (JijV j (t ~ njM ~ tw))) dt 
< dw l t =a a (b a v l t -w l t )dt (14) 
dyl = (A a S a (vl)(l-yl)-D a yl)dt + a Y (v,y)dBl> Y 

A similarly (but more complex) expression is obtained for the HH model using equations (12) 
and with distributed delays. 
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A. 2. Stochastic firing-rates models 

A phenomenological neuron model consists in considering that neurons interact through their 
mean firing-rate. The firing-rate model considers that the membrane potential has a linear 
dynamics, and its mean-firing rate is a smooth sigmoidal transform of the membrane potential 
S(r a , •) depending on the neural population a. In other words, an incoming firing rate provokes 
postsynaptic potentials that linearly sum. The neurons receive additional inputs that are the 
sum of a deterministic current I(r a ,i) and noise o~(r a )dWl . The network equations hence read: 

dV\t) = (-J—V\t)+I(r a ,t) + J2j^ E S(r 7 ,V j (t-T aj ))\ dt + <j(r a )dW l t 

\ [ra) 7=1 7 j,p(j)=1 ' 

It is easy to check that assumptions (H1)-(H5) are satisfied for the firing-rate model. 
Appendix B: Generalized models 

In the main section we chose to concentrate on the cornerstone mathematical problems arising 
in the modeling of neural fields, and chose to deal with relatively general models, yet simplified. 
Indeed, as discussed in section A, we saw that two technicalities were not taken into account in 
our general analysis. These were (i) non-globally Lipschitz drift that do not satisfy the linear 
growth condition (for Fitzhugh-Nagumo models) and (ii) state-dependent diffusion coefficients. 
Relatively classical methods allow to extend our proofs to these two refinements. In this section 
we will explain how one can extend to present analysis to models including this kind of dynamics. 
We consider the network equations: 

1 P{N) 1 

d Xi' N =(f(r a ,t,Xi' N )+" Yl E ^r^Xl'^X^^dt 

1 ' 7=1 p(j)=7 7 

+ 9{r a ,t,Xl' N )dWi. (15) 

hence make the following generalized assumptions: 

(HI') f(r, t, ■) and g(r, ■) are locally Lipschitz-continuous 

(H4') x T f(r,t,x) + \\g{r,t,x)\ 2 < K (1 + \x\ 2 ) uniformly in (r,t) 

The main difficulty is the non-global Lipchitz continuity of the drift and diffusion functions. 
However, under assumption (H4'), we can show that any possible solution is in M 2 ([— t, T], lL|(r)), 
thanks to the following lemma: 

Lemma 5. Let (Ct( r )) € A4 2 ([— r, 0], IL^(r)) be an initial condition for the mean-field equa- 
tion (3). Any possible solution (X t )t^[- T ,T] of the equation (3) with initial condition £° and 
measurable law with respect to B(T) is square integrable, in the sense: 

sup E [£ r \X t (r)\ 2 ] < C(T) (16) 

t&[-r,T] 

where C(T) is a quantity depending on the horizon T and the parameters of the system. 
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Proof. The proof is based on the application of Ito's formula for the squared modulus of X tl 
standard inequalities and Gronwall's lemma. In details, let X be a solution of the mean-field 
equations, and 

t„ = inf {t > ; ll^tHjL^r) > n). 

Due to the non-standard nature of the equation, let us underline the fact that Ito's formula is 
valid, i.e. that for any r G T, t X t (r) is a semimartingale. By definition, it is clear that for 
any r G T both X t+S (r) and Z t + S (r) are JF t measurable for all s G [— r, 0] since these are driven 
by a standard Brownian motion W t (r), implying that X t (r) is the sum of a continuous adapted 
process of finite variation: 

J (/(r, s, X s (r)) + £ r , [Mz[b(r, r', X s (r), Z s _ T(rif ,) (/))]) ds, 

a continuous (Tt, P)-local martingale which is a stochastic integral of a progressively measurable 
processes with respect to a Brownian motion, J g(r 7 s, Xf) dW s (r). 

We can therefore apply Ito's formula ( see e.g. [8, Chap. 4.5] in the context of delayed 
equations) to IE [|Xt Ar?i | 2 ], and obtain: 



£ r [Ep S Ar»| 2 ]] =E 



£ r [ |CoWr+2 sup 

se[o,t] Jo 



du { XZ(r)f(r,u,X u (r)) + -\g(r,u, X u {r))\ 2 



+ £ r , X^(r)^ z [b(r, r', X u (r),Z u _ T ^ rl) {r'))} 



<E[f r [|C (r)| 2 ]]+2(^ 



duM[£ r (l + \X u (r)\ 2 )] , 

' J o 

yielding, using Gronwall's lemma^ 

sup E [£[\X tAT J 2 ]] < E [£r[\C° (r)\ 2 ]] e K ' T =: C'(T) 
te[o,T] 

with K' = 2(K + \/Kk). This estimate is valid for any n, thus for n sufficiently large, the 
probability of \X t \ to exceed n prior to time T vanishes. Letting n go to infinity, we obtain: 



sup E 

te[-T,T\ 



\ X t\U : 



2 



< max sup E 

\s6[-t,0] 



,C'(T) =:C(T). 



□ 



The same result can be shown when considering the network equations. These results being 
proved, it is then possible to prove analogous versions of theorems 2 and 3. The state dependent 
diffusion function is easily controlled using Burkholder-Davis-Gundy theorem, and the proofs of 
theorems 2 and 3 are hence valid when considering the truncated drift and diffusion functions: 



and 



fu(r,t,x) 



gu{r,t,x) 



f(r,t,x) \x\<U 

J(r,t,Ux/\x\) \x\>U 

r g(r,t,x) \x\<U 

g(r,t,Ux/\x\) \x\>U 
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which are both globally Lipschitz-continuous, and from these results and lemma 5, it is folklore 
to extend these results to the original problem. 

Let us for instance focus on the existence and uniqueness of solutions for / and g functions 
that are not globally Lipschitz. Denoting Xjj the unique solution to the truncated problem, 
and defining the stopping time tjj — inf{£ 6 [0, T], ||A;y(i)|L 2 , r . > U}, it is easy to show that 

Xu(t)=Xu>(t) if < t < T V , U' > U, 

implying that the sequence of stopping times tjj is increasing. Using proposition 5 which implies 
that the solution to (3) is almost surely bounded, ensures existence and uniqueness of solutions 
of the mean-field equations (3) by letting U go to infinity. A similar argument applies for the 
propagation of chaos property. 
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